Raman signature of the U(l) Dirac spin-liquid state in spin-1/2 kagome system 
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We followed the Shastry-Shraiman formulation of Raman scattering in Hubbard systems and 
considered the Raman intensity profile in the spin-1/2 "perfect" kagome lattice herbertsmithitc 
ZnCu3(OH)gCl2, assuming the ground state is well-described by the U(l) Dirac spin-liquid state. 
In the derivation of the Raman T-matrix, we found that the spin-chirality term appears in the A^ g 
channel in the kagome lattice at the t 4 /(u!i — U) 3 order, but (contrary to the claims by Shastry and 
Shraiman) vanishes in the square lattice to that order. In the ensuing calculations on the spin-1/2 
kagome lattice, we found that the Raman intensity profile in the E g channel is invariant under an 
arbitrary rotation in the kagome plane, and that in all {A\ g , E g , and A2 g ) symmetry channels the 
Raman intensity profile contains broad continua that display power-law behaviors at low energy, 
with exponent approximately equal to 1 in the A-z g channel and exponent approximately equal to 3 
in the E g and the Ai g channels. For the A2 g channel, the Raman profile also contains a characteristic 
1/uj singularity, which arose in our model from an excitation of the emergent U(l) gauge field. 
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I. INTRODUCTION 

Quantum spin liquids, which are quantum ground 
states of interacting spin systems that break no sym- 
metries, have long been a fascination among the theo- 
retical condensed matter community. After years of ex- 
perimental searches, several promising candidates finally 
emerged, including the "perfect" spin-1/2 kagome lattice 
herbertsmithitc ZnCu3(OH)eCl2, which shows no signs 
of magnetic ordering down to a temperature of 50 mK, 
despite having a nearest-neighbor antifcrromagnetic ex- 
change J « 190 

Ever since the successful synthesis of herbertsmithitc, 4 
a host of experimental techniques have been ap- 
plied to study the material, including thermodynamic 
measurements neutron diffractionji 7 NMRpiL^ 
and /zSRi^ Unfortunately, the experimental results ac- 
cumulated thus far are still insufficient to determine if the 
material is truly a quantum spin liquid. In particular, the 
valence bond solid (VBS) state proposed in Ref. [9land[l0l 
remains a possible alternative to the U(l) Dirac spin- 
liquid (DSL) state proposed in Ref. [ill In order to set- 
tle the debate of which theoretical model best describe 
the quantum state in herbertsmithite, further experimen- 
tal probes, with guidances from theory, are probably re- 
quired. 

Recently, Cepas et alX^ considered Raman scatter- 
ing on the spin-1/2 kagome system and concluded that 
a generic spin-liquid state can be distinguished from a 
generic VBS state by the polarization dependence of the 
signal. They also obtained a more detailed prediction of 
the Raman intensity using a random phase approxima- 
tion, which may be too crude given the subtle orders 1 ^ 
that may be present in the system. 

In this paper, we consider the specific scenario that 
herbertsmithite is described by the DSL state, and make 
theoretical predictions of the experimental signatures 
that should be present in the Raman scattering. In con- 
trast to the predictions by Cepas el al. (who argued that 



signals consist of several sharp Lorcntzian peaks), wc pre- 
dict that broad continua should be present in all (E g , A\ g , 
and A2 g ) symmetry channels, with a super-linear power 
law at low energy in the A\ g and the E g channels, and 
a linear power law in at low energy in the Ai g channel. 
In addition, we predict that ol/u singularity should be 
observed in the Ai g channel. From a theoretical perspec- 
tive, this 1/u) singularity is particularly interesting, since 
it can be thought of as arising from an excitation associ- 
ated with an emergent gauge boson in the system. 

The paper is organized as follows: In Sec. [H] the 
Shastry-Shraiman formulation of Raman scattering in 
Hubbard systems is reviewed and our disagreement with 
the original results by Shastry and Shraiman on the 
square lattice is discussed, after which the relevant re- 
sults for the kagome lattice is presented (the detailed 
derivations and further discussions are relegated to Ap- 
pendices). In Sec. PH the U(l) Dirac spin-liquid model 
is reviewed, with important features in the mean-field 
theory highlighted. The Raman intensity profile in E g 
channel is then presented in Sec. IIV| and the analogous 
results for the A\ g and A2 g channel are presented in 
Secs.lVland lVIl respectively Further discussions on these 
results are presented in Sec. IVIII 



II. SHASTRY-SHRAIMAN FORMULATION 

Being a strongly correlated material, herbertsmithite 
can be described using a one-band Hubbard model: 

#Hb = H t + H v = - Uj4a c j<y + U ^2 n>it n H i (!) 



where i, j label lattice sites and a =t, 4- labels spin, cj 
(cj) is the electron creation (annihilation) operator on 
site i, and n i(7 = c\ a c la . 

Coupling to the external electromagnetic field 
can be incorporated by the replacement c\ a c, J(y i— > 
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q a .Cj (T exp(|| J? A • dx). Expanding this exponential 
and including also the free photon Hamiltonian H-y, the 
Hamiltonian now reads: 

H = i?Hb + H 1 + He , 

e 2 / .x,+x,, . ,\2 \ ( 2 ) 

q 

where (a^) denotes the photon creation (annihila- 
tion) operator at momentum q and polarization a, and 
A(x) denotes the photon operator in real space. The • • • 
are terms at higher order in A. 

By treating He as a time-dependent perturbation, the 
transition rate from an initial state |i) to a final state |/) 
is given by: 

T fl = 2*\{f\T\i)\ 2 8{S f -S, l ) , (3) 

where Si (£/) is the energy of the initial (final) state and 
T = H c + H c (£i - H Hh -H 7 + it))- 1 H c +■•• is the 
T-matrix. 

Since the fine-structure constant e 2 /hc « 1/137 is 
small and since we are interested in Raman processes 
(one photon in, one photon out), only terms second or- 
der in A need to be retained. At this order, the T-matrix 
reads: 

T = H { J } +H^ £i _ {H J + ^ ) + i7] H^ = T NR + T R , 

(4) 



Next, a spin quantization axis is fixed and the initial 
states |i) = \{<r}) <8 |k»,ej) and final states \ f) = \{cr'}) ® 
|kj, e/) are taken to be a direct product of a definite spin 
state in position basis with a photon energy eigenstatei 2 ^ 
Then, a complete set of states is inserted in between the 
operators in Eq. [5] By the assumption U 3> U\ > \t\, 
the intermediate states are dominated by those having 
no photons and exactly one holon and one doublon. 
Thus, they take the generic form jr^; r^; {r}) CS> |0) , where 
\rd;r h ;{r}) = 4 dl< 7 c n.,o-)l{ T }) is obtained from the 
spin state \{t}) by removing an electron at and 
putting it at r^, and |0) denotes the photon vacuum state. 
Henceforth we shall take the abbreviation that spins are 



where H c denotes the part of He that is n-th order in 
A. The subscript R and NR on the last equality stands 
for resonant and non-resonant, respectively. 

We are interested in a half-filled system ((J2a n i<?) = 
1) in the localized regime (U ^> i), in which both the 
initial and the final state belongs to the near-degenerate 
ground-state manifold n^n^ = 0. In such case, Tnr has 
no matrix element that directly connects between the 
initial and the final states. Hence, only Tr is relevant for 
our purpose. 

Let u>i (idf), ki (k/), and (e/) be the frequency, 
momentum, and polarization of the incoming (outgo- 
ing) photon, respectively. Then, Si = + S^' = 
LJi + 0(t 2 /U), where H Hb \i) = s!> Hh) \i); and A(x) i-> 
g^iCL^.e^ * + g f efa^e'' k f x , where gt = ^hc 2 /io kl fl 

and gf = yjhc 2 /w^fl, with VL being the appropriate vol- 
ume determined by the size of the sample and / or the size 
of the laser spot. In much of the following we shall assume 
as typical that the momenta carried by the photons are 
much smaller than the inverse lattice spacing, and hence 
g-tki-x ^ g-ik/ x ^ ghgjj a i so assume that the sys- 
tem is near resonance, so that U 3> |o;, — U\ > \t\. Con- 
sequently, henceforth we shall keep only terms that are 
zeroth order in t/U and expand in powers of t/(wi — U). 

Since the initial and final states both belong to the 
near-degenerate ground-state manifold, it should be pos- 
sible to re-express Tr in terms of spin operators. A 
procedure for doing so was developed by Shastry and 
Shraimani i 15 A first step in the derivation is to expand 
the denominator of T R : 



(5) 



I 

summed implicitly within pairs of electron operators en- 
closed by parentheses, so that, e.g., (c\cj) = ^ a c*i a Cja- 

Under this insertion, (Si — H;j — H^)^ 1 = (a;, — U)^ 1 
becomes a c-numbcr. Moreover, recall that H t and (ne- 
glecting the photon part) He arc sums of operators of the 
form (c\cj). Once a particular term is picked for each of 
these sums, and given an initial spin state |{cr}), the re- 
sulting chain of operators automatically and uniquely de- 
termines the intermediate states (which may be 0). Thus 
the intermediate states can be trivially re-summed, and 
Eq. O becomes, in schematic form: 
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The sum in Eq. [6] is formidable. However, if He and 
H t connects only between sites that are a few lattice con- 
stants away, then at low order in t/(uji — U), except for 
the choice of the initial site (ji in Eq. the number of 
non-zero terms is finite and does not scale with the lattice 
size. Thus, Eq. [6] provides a systematic way of analyzing 
the contributions to the Raman intensity. 

The final step in the Shastry-Shraiman formula- 
tion is to convert the chain of electron operators 
( C I c j n )'" (cjj c ji ) m to spin operators using the anti- 
commutation relation and the following spin identities: 



1 



f - y , = -5 , 



(7) 



S- 



where S = cj. (t ctct < /2)c ct < is the spin operator for spin-1/2 
and t is the usual Pauli matrices. 

To the lowest non- vanishing order in t/(uii — U), the 
Shastry-Shraiman formulation reproduces the Fleury- 
London HamiltonianJ^ i.e.: 
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FIG. 1: Pathways that contribute to the spin-chirality term in 
the square lattice. A thick (blue) arrow indicates the initial 
or the final hop, which arises from He- A thin (magenta) 
unbroken arrow indicates the movement of a doublon, and 
a thin (magenta) broken arrow indicates the movement of a 
holon, both arising from Ht. 



(f\T R \i) = {{a'}\H F ^\{a}) + 
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H FL = T r (et • M)(e/ ■ M 
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U — u>i 



(8) 

where [i = r' — r is the vector that connects lattice site 
r to lattice site r'. 

For the square lattice with only nearest-neighbor hop- 
ping, Shastry and Shraiman claimed that spin-chirality 
term Si ■ (Sj x Sfe) appears in the ejef — e^ef channel 
at the t /{iOi — U) 3 order. However, our re-derivation 
does not confirm this result and instead concludes that 
the spin-chirality term vanishes to this order. For details 
about our re-derivation, see Appcndix[A] While it is hard 
to pin down the source of this discrepancy, two possibili- 
ties are plausible. First, the pathways that contribute to 
the spin-chirality term includes not only those in which 



a doublon or holon hops through a loop [Figs. 1(a) and 
l(b)| , but also those in which a holon "chases" a doublon 
or vice versa without involving a fourth site [Figs. 1(c) 
and l(d)| . These chasing pathways are non- intuitive and 



could be easily missed. Second, observe that for c\c, 
the spin indices are flipped in Eq.[7|when going from elec- 
tron operators to spin operators. This, together with the 
applications of the anti-commutation relation, can easily 
produce minus sign errors. 

Our conclusion that the spin-chirality term vanishes to 
the t A I (u>i — U) 3 order in the one-band Hubbard model 



need not contradict with the experimental claim that the 
spin-chirality term has been observed in the cuprates^ 
for in the cupratcs — with the holon being delocalized as 
Zhang-Rice singlet while the doublon being localized at 
the copper site — the holon and doublon hopping magni- 
tude need not be equal. In that case, the crucial can- 
cellation between the four pathways in Fig. Q] no longer 
occurs. Furthermore, in the Shastry-Shraiman formal- 
ism the spin-chirality term may also be present at higher 
order in t/(u)i — U) and/or when further neighbor hop- 
pings are included. Since the ratio t/(u)i — U) need not 
be small near resonance, these higher-order effects can 
manifest in experiments. 

Now we specialize to the kagome lattice and for sim- 
plicity assume that the hopping is between nearest neigh- 
bors only. Our convention of lattice basis, primitive lat- 
tice vectors, and axis alignment is shown in Fig. [2] In 
contrast to the square lattice, we found a non-vanishing 
spin-chirality term at the t 4 /(uji — U) 3 order. The con- 
trasting result between the kagome lattice and the square 
lattice can be traced back to the lack of four-site loop in 
the former. 

For theoretical calculations, it is convenient to decom- 
pose the polarization dependence of the Raman inten- 
sity into the irreducible representations (irreps) of the 
lattice point group, since operators belonging to differ- 
ent irreps do not interfere with each other (note however 
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FIG. 2: The kagome lattice in real space, with the unit vectors 
x, y, primitive lattice vectors ai, &2, and lattice basis 1, 2, 3 
defined. 

that subtractions between various experimental setups 
are often required to extract the signal that corresponds 
to a particular channel^). It is knowr>i& that herbcrt- 
smithitc belongs to the space group R3m and hence to 
the point group In D^d, the polarization tensor 

S Q B=x y ^a,/3^/ e f m l ne kagome plane decomposes into 
two one-dimensional irreps A\ g and A2 g , and one two- 



dimensional irrep E g : 
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+ e v f e x 



(9) 



To the lowest non- vanishing order in t/(u>i — U), the 
Ai g and the E g component of the T-matrix are derived 
from the Flcury-London Hamiltonian Eq. [5J However, 
the resulting expression for the A\ g channel is the sum 
of a constant and a term proportional to the Heisenberg 
Hamiltonian and thus at zeroth order in t/U does not 
induce any inelastic transitions. For inelastic transitions 
in the Ai g channel, the leading-order contribution ap- 
pears at the t 4 /(uJi — U) 3 order instead. The leading- 
order contribution to the A2 g channel also appears at the 
t 4 /(uji — U) 3 order. Explicitly, to leading order and ne- 
glecting the elastic part, the operators that corresponds 
to the different channels are: 



Oa 2 



4f- 



w< - U ^ V 4 



( tS_R,3 • (Sr,1 + Sii +a2 .i + S fl: 2 + Sfl+aa-ai.a) - ^Sj?,2 1 (Sfl,l + Sfl+ai.l) J , 



4/- 



uji-U ^ 4 V 

R 



(10) 

(11) 



— t ( 

OAi s = -j— ^3 2 ( Sr,1 ■ (SiJ+ ai ,l + S j R +a2 ^i) + ■ (S_R_ ail 2 + S^ +a2 _ ai . 2 ) + Sr^ ■ (S#_ a2 .3 + Sfi_ a2+ai ^ 3 ) 



Sfl,l • (Sf{_|_ a2 _ aij 2 + Sfl +ai _ a2r 3) + Sft t 2 ■ (Sfl-a 2 ,3 + Sfl+ a2i i) + Sfl,3 ■ (Sfl+ ai ,l + Sfl_ ai ^)) , (12) 

2\/3it 4 



(ui - uy 



^ (3Sr,1;R,2;R,3 + 3<SH,l;fi-ai,2;R-a 2 ,3 + <->R,l;-R,3;B+a 2 -ai,2 + ■SR-faa,!;^;^ 



Sr,3 ;i?,2;i?+ai ,1 + 5ji_ a2 +a 1 ,3;ii 1 2;fl,l + Sr,2;R,1;R-* 2 ,3 + <Sfl-a l ,2;i?,l;i? v 



{uji - uy 



?( 3 A +3 V+t 



(13) 



where Si-j-k denotes Si • (Sj x Sfe). Further discussion on 
the derivation of these expressions can be found in the 
Appendices. 



the Hubbard model is reduced to the Heisenberg model: 



1 Si • Sj 



1 



(14) 



III. U(l) DIRAC SPIN-LIQUID STATE 



At half filling 



1), on the near-degenerate 



ground-state manifold and to the leading order in t/U, 



where Jy = Atf^/U . As in the previous section, we shall 
specialize to the case where the hopping is between near- 
est neighbor only, and hence Jy = unless i, j are near- 
est neighbor. 

In a spin-liquid state (Si) = 0. Thus, new variables are 
needed to describe the order of the system. A common 
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Energy (units of \J) 




(a) 



(b) 



FIG. 3: (a) The kagome lattice with the DSL ansatz. The 
(red) dashed lines correspond to bonds with effective hopping 
t = — xJ while the (blue) unbroken lines correspond to bonds 
with effective hopping t = \J . a' x and SL2 are the primitive 
vectors of the doubled unit cell, (b) The original Brillouin 
zone (bounded by unbroken lines) and the reduced Brillouin 
zone (bounded by broken lines) of the DSL ansatz. The dots 
indicate locations of the Dirac nodes at half-filling, the crosses 
indicate locations of the Dirac nodes crossing the second and 
the third band, and the thin (gray) lines indicate saddle re- 
gions at which the band energies are the same as that at 
k = 0. 

choice is to introduce fermionic spinon operators , fj a 
to decompose the spin operator Si as follows: 



t T 



pfi/3 



(15) 



The occupation constraint rii a = 1 resulting from the 

U — > oo limit then becomes the constraint f} a fia = 1 
on the spinonsi 19 ' 20 

Under this transformation, the spin-spin interaction 
Si ■ Sj is mapped to a four-fermion interaction, which 
can be decoupled using a Hubbard-Stratonovich trans- 
formation. This yields the partition function Z = 

J DfDpD\D X cxp (- Jjf drL\ , with the La grangian L 

given by: 



{ij) \ <y 



(16) 



Mean-field ansatzs can be specified by treating Xij as 
an order parameter. Since a spin-liquid state is invariant 
under translation and rotation, \xij \ must be independent 
of i, j. Letting Xij = X e ~ laij an d rewriting Aj = Oq, the 
Lagrangian Eq. [TH] yields the following mean- field Hamil- 
tonian: 

ia (ij),t? 

(17) 
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FIG. 4: The band structure of the kagome lattice with the 
DSL ansatz plotted along (a) k x = and (b) k y — 0. Note 
that the top band is twofold degenerate. 



Observe that an internal gauge field emerges nat- 
urally from this formulation. Its space components ay- 
arise from the phases of Xij > while its time component ao 
arises from enforcing the occupation constraint. 

By gauge invariance, a mean-field ansatz for a M is 
uniquely specified by the amount of fluxes through 
the triangles and the hexagons of the kagome lat- 
tice. In particular, the DSL state is characterized by 
zero flux through the triangles and n flux through the 
hexagons . 9 i 11 ' 13 By picking an appropriate gauge, the 
DSL state can be described by a tight-binding Hamil- 
tonian with effective hopping t = ±x</ on each bond. 
For the precise pattern of the ± signs see Fig. 3(a) Note 
that the unit cell for this effective tight-binding Hamil- 
tonian is necessarily doubled as the flux enclosed by the 
original unit cell is n. 

It is easy to check that, in units where xJ = lj this ef- 
fective tight-binding Hamiltonian produces the following 
bands: 



-Stop — 2 



(doubly degenerate) 



(18) 



E 



1± a/3=f V2\/3-cos2fc 



2 cos k x cos Vofcj, . 

(19) 



At any k-point, < £L _ < E + _ < E + _ + < E top . 



For plots of this band structure, see Fig. [4] 

At low energy, the spinon spectrum is well-described 
by four (two spins times two k-points) Dirac nodes, lo- 
cated at momentum ±Q = ±n/\/3y [c.f. Figs. |3(b~)1 and 
4(a)| . More specifically, at low-energy we may replace the 
mean- field Hamiltonian Eq.[l7]by the Dirac Hamiltonian: 



#Dirac = V F ^ i } l,a,q{ ( lxTx + q y Ty)lpa, a/ 



(20) 



where a =t, i index spins, a = ± index the location of 
the Dirac node, and q denotes the momentum as mea- 
sured from the Dirac node. The relation between the 
two-component fermionic operators tp, %jy and the spinon 
operators f,p defined on the lattice sites can be found 
in Ref . [H 
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As will be explained in details below, the Dirac struc- 
ture of this low-energy Hamiltonian bears important con- 
sequences for the Raman intensity at low energy. In par- 
ticular, the power-law behavior of the Raman intensity at 
low energy is largely resulted from the strong phase space 
restriction of the Dirac node structure. 

Under the DSL ansatz, we may take \i) = |^ Hb ^) ® 
|kj,ej) and |/) = |/' Hb ^) <B> |k/,e/) in the transition rate 
Eq. [3j where |i( Hb )) and |/ ( - Hb ' 1 ) are states obtained from 
filling the spinon bands. In particular, at zero tempera- 
ture (which will be assumed henceforth), |i( Hb )) is simply 
a state with spinons filled up to the Dirac nodes at the 
top of and |/^ Hb ') are states with a few spinon- 

antispinon pairs excited from |«( Hb )). 

Moreover, the spin operators that appear in Eqs. [TOL 
1131 can be converted to spinon operators using Eq. 1151 
Explicitly, 



Sj • Sj — ^ fiafj&fja'fiu' j 



(21) 



Si ■ (S, x S fc ) = \ {[flf 3 ){f}h)(flf l ) - h.c.) , (22) 

which allows the matrix elements between |^ Hb ^) and 
|/( Hb )) to be calculated at mean field using Wick's theo- 
rem. 

At zero temperature, to obtain the overall Raman in- 
tensity I a for channel a at a given Raman shift Aw = 
LUi—ujf, all final states that satisfies the energy constraint 



£ 



(Hb) 
/ 



£ 



(Hb) 



Aw must be summed over. Strictly en- 



forcing this constraint is difficult in a numerical computa- 
tion. Instead, we sample the final states |/ ( - Hb ' ) ) without 
imposing the energy constraint, but sort them into bins 
of energy. The overall intensity is then obtained by sum- 
ming all states whose energy fall within the same bin. 
We perform the summation numerically using a simple 
Monte Carlo sampling of the momenta and band indices 
of the spinon and antispinon excitations in |/^ Hb ^). 

Some calculations are also performed in a slightly dif- 
ferent way, by first converting the summed squared am- 
plitude into a correlation function: 



\(fW\O a \i 



(Hb) \ 



'5{£f h) 



r(Hbj \ \ 

c i ~ L ^ ul ! final states explicitly Given the form of 0„(i) and O , 



= / dte iAut {i^\O a (t)O a (0)\i^) 



/a (Aw) 



(23) 

To calculate this correlation function, the spinon opera- 
tors that appears in O a (c.f. Eqs. [JU][T3j are converted 
from real space to momentum space, which introduces 
sums over momenta and band indices. Such sums are 
again computed by simple Monte Carlo samplings in the 
manner explained above. 

To simplify notations, in the following we shall abuse 
notation and write \i) in place of |^ Hb ^) and similarly 
write |/) for |/^ Hb ^) whenever the context is clear. 



IV. E a CHANNEL 



First consider the Raman intensity in the E g channel, 
Ie = 1^(1) + ■ By computing the correlation func- 

tion Eq. [23] (with a = Eg 1 ^ and Eg), we obtain the 
Raman intensity profile as shown in Fig. [5] 

From Fig. [5] it can be seen that the Raman response 
take the form of a broad continuum ranging from approx- 
imately 1.5x^ to approximately WxJ , with occasional 
sharp spikes. The existence of a continuum is not a sur- 
prise, since the operators 0„{i) and 0„{2) in Eqs. 1101 - 
1111 corresponds to two-spinon-two-antispinon operators 
in the DSL ansatz, and hence the final states |/) arc 
spinon- antispinon pairs and thus there is a continuum of 
phase space for excitations. The cutoff near WxJ i s a l so 
natural, since the total band-width in the DSL ansatz 
is approximately 6 \ Jy an d hence with two spinons and 
two antispinons the excitation energy is at most 12 \J- 
The sharp peaks and the low-energy suppression, how- 
ever, require further investigations. 



Intensity (arb. units) 




4 6 8 

Energy shift (x J) 

FIG. 5: Raman intensity in the E g channel, computed using 
the correlation Eq. 1231 with a bin size of O.Q5\J m ener gy- 



To gain more insight into how the various features in 
Fig. E3 comes about, it is useful to consider the sum over 



in Eqs. [10] [Tl] there are two types of final states: the first 
consists of two spinon-antispinon pairs excited from the 
ground state \i) while the second consists of one spinon- 
antispinon pair excited from \i). Schematically, these 
take the form: 

|/2pairs) = /fc ll n 1 ,CT^2,«2,ff/fc 3l n3,a-'/' i: 4>™4> '' I*) ' (^) 
I /l pair) = fk ltni .afk2,n 2 -<y\i) ■ (25) 

In both cases the momentum conservation odd kj = 

Ek, holds, and that the band index n,- denotes one 
i even L L 

of the top three (empty) bands when i is odd and one of 
the bottom three (occupied) bands when i is even. 
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Intensity/DOS (arb. units) 



> i 




2 4 6 8 10 

Energy shift (xJ) 



FIG. 6: Plots of density of states (thin red curve) and E g 
Raman intensity (thick blue curve) contributed by final states 
having two spinon-antispinon pairs, together with the overall 
E g Raman intensity (broken gray curve), all computed with 
a bin size of 0.05xJ in energy. Note that the relative scale 
between the intensity plots and the DOS plot is arbitrary and 
is set here such that both curves are visible. 
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1 2 3 4 5 6 
Energy shift (x-0 



FIG. 7: Plots of density of states (thin red curve) and E g 
Raman intensity (thick blue curve) contributed by final states 
having one spinon-antispinon pair, together with the overall 
E g Raman intensity (broken gray curve), all computed with 
a bin size of 0.05x»7 in energy. Note that the relative scale 
between the intensity plots and the DOS plot is arbitrary. 



For the one-pair final state to have a non-zero matrix 
element with the initial state, a pair of spinon-antispinon 
operator must self-contract within 0„(i) and O p (2>, e.g., 

9 9 

fiafjvfja' fin' ^ (fiafj^fja'fiv' + fia fj" (fjcr' fif') 

~ (fiafi<?')fja'fj<? ~ (fja'fj^fiafic' 

~ Xijfjcrficr + Xjifiofjo ~ ^{fiofic "h fjafjc) ' 

(26) 

where the spin indices a, a' are summed in the above. 
The mean-field parameter Xij here is the same as the 
one introduced in Sec. Mil 

However, one should be mindful that while the con- 
straint fjpfjcr — 1 is strictly enforced in the exact 
spin-liquid state, it is enforced only on average in the 
mean-field representation. Therefore, while operators of 
the form fj a fj a should induce no transition between \i) 

and \f), the matrix element (f\fj a fja\i) may be non-zero 
at mean-field. To avoid such problem, we shall throw 
away the fj a fja terms in Eq. [26] by hand. 

With this precaution, we recalculated Ie 9 by numeri- 
cally sampling the final states, and we do so for the one- 
pair and two-pair contributions separately. The results 
are shown on Figs. [6] and [7] We have verified that the 
sum of Raman intensities from the two figures (with the 
relative ratio determined by the details of the numerical 
calculations) produces an overall intensity profile that 
matches Fig. [5] 

The density of states (DOS) for two-pair and one-pair 
excitations having zero total momenta have been plot- 
ted alongside with the respective Raman intensities in 



Figs. [5] and [7J From the figures, it can be seen that 
the DOS matches the Raman intensity profile very well 
for two-pair excitations, and less so (but still reason- 
ably well) for one-pair excitations. This can be under- 
stood by rewriting the first line of Eq. [53] as 7 Q (Aw) = 

\/|O a |£) T>(Aw), in which T>(£) denotes the density of 

state at energy £ and \/|O a |i) denotes the average 
matrix element squared at the same energy. From this, 
the DOS is expected to match the Raman intensity well 
as long as the average matrix element does not change 
drastically with energy — an assumption more valid for 
two-pair states as opposed to one-pair ones, because of 
the larger phase space available in the former case. 

Moreover, the sharp peak appearing near 5.5x>/ in 
Fig. [5] can now be attributed to one-pair excitations, 
which can in turn be attributed to a peak in the one- 
pair DOS. From the band structure (Fig. [4j, it can be 
checked that Aw = 5.41% J corresponds to the energy dif- 
ference between the top flat band and the saddle point 
k = at the bottom band. The enhanced phase space 
near this energy is thus the likely cause of this peak-like 
feature. 

Next, we consider the low-energy (Aw < 2\J) part 
the of Raman intensity more carefully, checking if the 
intensity profile is an exponential or a power law, and 
determining the exponent if the latter case holds. To 
enhance the quality of the data, and since only the low- 
energy behavior concerns us, we resample the final states, 
restricting the antispinon to the highest occupied band 
E + _ and the spinon to the lowest unoccupied band E +t+ 
(despite this, the two-pair data is consistently non-zero 
only above Aw = 0.4% J). The resulting data arc shown 
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2Q-0.25 -lqO.OO 10 - 25 

Energy shift (xJ) 

FIG. 8: Log-log plot of the DOS (red x symbols) and Raman 
intensity in the E g channel (blue + symbols) for two-pair 
excitations. Simple linear fits (straight lines in pink and cyan) 
give a slope of 5.1 for the DOS data and 4.9 for the Raman 
intensity data. 




Energy shift (xJ) 

FIG. 9: Log-log plot of the DOS (red x symbols) and Raman 
intensity in the E g channel (blue + symbols) for one-pair 
excitations. Simple linear fits (straight lines in pink and cyan) 
give a slope of 1.0 for the DOS data and 2.8 for the Raman 
intensity data. 



in the log-log plots in Figs. [8] and OH Since the smaller 
exponent dominates as Aw — > 0, the overall Raman in- 
tensity in the E g channel scales roughly as oc (Aw) 3 
at low energy. 

From the plots, it is clear that the DOS and the Raman 
intensity both follow a power law, with a higher exponent 
for the two-pair excitations as compared to the one-pair 
ones. Analytically, it is easy to check that 2?i pa ir oc Aw 
for the Dirac Hamiltonian Eq. since the integral in- 
volves four components of momenta subjected to three 
constraints (energy and momentum conservation), and 
that |q| scales as energy. Similarly, it is easy to check that 
2?2pair oc (Aw) 5 in the Dirac Hamiltonian. Furthermore, 
since the eigenstates of the Dirac Hamiltonian depend 
only on 9 = tan -1 (q y /q x ) but not on |q|, the average ma- 



trix element squared \ (f\O a \i)\ must be constant in en- 
ergy. However, to the t 2 / (w^ — U) order in the E g channel 
(Eqs. fT0l - [TT|l . the matrix element turns out to be exactly 
zero for all one-pair excitations in the Dirac Hamiltonian. 
Hence, at low energy, we expect -?2pair oc Impair oc (Aw) 5 
and /i pa i r oc (Aw) Q with a > 1, consistent with the nu- 
merical results in Figs. [5] and [5J 

Since 2?i pa ir oc Aw, it leaves open the possibility that 
Ie oc Aw at higher order in i/(wj — U). While we are 
not able to rule out this possibility, we do find that the 
vanishing of matrix elements in the Dirac Hamiltonian 
for all one-pair excitations persist to the t^jiioi — U) 
order (the derivation of Oe 9 to this order can be found 
in Appendix [Bjj^i Hence, even if Ie oc Aw at some 
higher order, its effect will not be prominent unless the 
system is sufficiently near resonance. 

In the data presented above we have summed the con- 
tributions arising from 0„(i) and 0„{2). By computing 

the two contributions separately, it can be checked that 
each contribute equally. In fact, we have checked that 
the intensity profiles are essentially identical upon an ar- 
bitrary rotation in the kagome plane for the one-pair and 
two-pair excitations separately. In other words, we found 
that the quantity /^(Aw,^), defined by: 

I Eg (Acj,e) = \{f\0 B n cos6 + Em sine\i)\ 2 

f (27) 
x«5(4 Hb) -ff b) -Aw), 

is independent of 9, and remain so even if the sum is 
restricted to one-pair or two-pair states. Our numerical 
results arc thus consistent with the analytical arguments 
given by Ccpas et al^ 

V. A lg CHANNEL 

Using Oai s m Eg. [T2] in place of O e (d and O e (2) , we 
repeat the calculation of the Raman intensity profile for 
the Aig channel. The results are shown in Fig. \W\ From 
the figure, we see that the Raman intensity profile of the 
Aig channel also has a broad continuum up to a cutoff 
near llxJ- However, the sharp peak near 5.5x</ that 
appears in the E g channel is markedly missing. 

Decomposing the Raman intensity into one-pair and 
two-pair contributions as in Sec. IIV1 it can be seen that 
the overall Raman intensity profile in the Ai g channel is 
dominated by the two-pair states. Since the sharp peak 
near 5.5%J is originated from the one-pair contribution, 
this explains the absence of sharp peak in the A\ g chan- 
nel. 

However, at low energy (Aw < 1.5% J) the Raman in- 
tensity profile is still dominated by the one-pair contri- 
bution. And by plotting the Raman intensity profile in a 
log- log scale (Fig. we see that the low-energy behav- 
ior is characterized by a power law with exponent a ~ 3, 
similar to the value obtained in the E g channel. Again, 
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Intensity (arb. units) 




Energy shift (xJ) 

FIG. 10: Plots of the overall Raman intensity (thick and blue) 
and its contribution by one-pair states (thin and red) in the 
Ai B channel, both with the same vertical scale and computed 
with a bin size of 0.05xJ in energy. 




Energy shift (%J) 

FIG. 11: Log-log plot of the Raman intensity (blue + sym- 
bols) and its contribution by one-pair states (red x symbols) 
in the Ai g channel. Simple linear fit (straight line in cyan) 
for the one-pair data up to Aoj = \ J gives a slope of 3.1. 



it can be checked on analytical ground that the ma- 
trix element vanishes for all one-pair state in the Dirac 
Hamiltonian, consistent with the numerical results. 



VI. A 2g CHANNEL 

For the A\ g channel, to the leading order in t/(oj% — 
U), the Raman process receives contributions from ex- 
cited states having one, two, and three pairs of spinon- 
antispinon. The overall Raman intensity coming from 
these spinon-antispinon pairs are plotted in Fig. 1121 It 
can be seen that the continuum also appears in this chan- 
nel, now ranging from 0%J up to approximately 16% J, 
which corresponds to approximately three times the total 
spinon bandwidth. Moreover, sharp peaks not unlike the 
one in the E g channel are observed at various energies. 



Again, it can be checked that these sharp peaks can be 
attributed to the various features of the DSL bands, par- 
ticularly to the flat-band-to-saddle and saddlc-to-saddlc 
transitions. 



Intensity (arb. units) 
> i 




2 4 6 8 10 12 14 16 
Energy shift (x J) 

FIG. 12: Raman intensity in the Ai g channel arising from 
spinon-antispinon pairs (thick and blue), and its contribution 
from two-pair (thin and red) and three-pair (thin and green) 
states. All computed with a bin size of 0.05x>/ in energy. 



Note that the low-energy Raman transition is much 
more prominent in the A 2g channel than that in the E g 
and A\ g channels. In particular, the Raman intensity 
in this channel has a broad peak near Aw = l.b\J , be- 
low which it is visibly linear. From Fig. 1121 the promi- 
nence of low-energy Raman transition trace back to the 
relatively large average matrix elements in the one-pair 
transitions. Moreover, the linear behavior suggests that 

\\f\O a \i)\ no longer vanishes in the Dirac Hamiltonian, 
so that Iipair oc 2?i P air oc Alo at low energy 

The forgoing discussion neglected an important con- 
tribution to the A2 g Raman intensity. Recall that in 
computing the one-pair contribution, a mean-field fac- 
torization is performed (Eg. I26[) . in which the Hubbard- 
Stratonovich variable Xij is treated as a constant. How- 
ever, because of the emergent gauge structure in the U(l) 
spin-liquid theory Eq. [17l the dynamics of the phase of 
Xij cannot really be neglected. Fortunately, since the 
projection of the T-matrix onto the E g and A\ g chan- 
nel Eqs. IT0HT21 involves no non-trivial closed path, the 
results presented above should still be qualitatively cor- 
rectly, even though quantitative corrections to the de- 
tailed predictions (such as the exponents of the power 
laws at low-energy) may be present. 

The situation for the A 2g channel is different, as can 
be seen by considering the full contraction of the spin- 
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chirality term Si • (S2 x S3) into Xij- 

iSr • (S 2 x S 3 ) = \ ({f\h){ftf2){flh) - h.c.) 

= (Ce 4(ai3+a32+Q2l) - h.c.) + . . . , 

(28) 

where Xy denotes the part of Xij that can be treated 

as constant, and C — Xi^X^X^iV^- For convenience we 
shall denote (C e i( ~ a ' k+a ^ +a ^ - h.c.) as Q ijk henceforth. 

Note that e^ ai3+a32+a21 ) is the emergent gauge flux en- 
closed by the loop l43-)24 1. Since the emergent 
gauge field is fluctuating, this gauge flux can lead to an 
excitation of the kagome system. Physically, the excita- 
tions generated by Q can be thought of as a collective 
excitation in the system, which exists on top of the indi- 
vidual spinon-antispinon excitations. This is analogous 
to the situation in an ordinary Fermi liquid, where plas- 
mon mode exists on top of the electron-hole continuum. 



Since the fluctuation of can be considered as a col- 
lective excitation in the system, as an approximation we 
may consider the final states |/ ga ugo) connected to ground 
state by Q to be lying in a separate sector than the 
spinon-antispinon pairs previously considered. Then, the 
total Raman intensity in the Ai g channel can be obtained 
by adding the contributions by these collective states to 
the contributions by the spinon-antispinon pair states. 

The contribution to the Raman intensity by |/ gaU ge) 
can be computed from the QQ correlator. To do so, we 
start with Eq. I28( take the continuum limit, and Tay- 
lor expand the exponential. Then, 013 + + «2i l— ^ 
J 132 a ' ^ X = -f-fiii32 ^ ^ 2x ' where- b = d x a y — d y a x 
is the emergent "magnetic" field. Here P132 denotes 
the closed loop 1 — > 3 — > 2 — > 1 and de- 
notes the area enclosed by this path. Thus, Qijk ~ 

2i\C\ (sinOo) + cos(>o) JJ Q _, k bcPxj w 2i|C|(sin(0 o ) + 

cos((/)o)Qijk b(rijk)), where C = \C\e l ^° and denotes 
the position of the three-site loop^ The QQ correlator 
is thus converted into a 66 correlator. Hence: 



7 ( ? u g e) ^ C ,J dteiAult{i{ lKrjilu .y (j2 n J2 A , s b(r ijk ,0)e-^) \i) , (29) 



where C is a numerical constant, and will not be kept 
track of below. Note that the photon momentum has 



been restored, which introduce the factor e 



'«* , with 



q = k; hf the momentum transferred to the lattice. 
The X^A's * s a shorthand for summing over the different 
three-site geometries on a unit cell with the appropriate 
coefficients, as is shown in Eq. [T3l 

Since 6 = d x a y — d y a x , 6(x) = J2k(^ x e fc) a ke lkx — 
(ik x ejl)a ] { ( .e _lk ' x , where is the polarization of the 
emergent gauge field at momentum k. Moreover, a,k\i) = 
and (i\a\. = 0. Hence, upon Fourier transform, 



/(f 2 a ; g0) « (i|(<q • e^AwJat (0)(-*q • e*)|i> 



(30) 



The correlator that we need to compute is thus one of 
a gauge field in 2+1 dimensions coupled to relativistic 
fermions (i.e., fermions described by the Dirac Hamilto- 
nian Eq. [20)) . This situation has been considered by Ioffe 
and Larkin^i under the context of high-T c superconduc- 
tivity, who found that: 

ll afi {q,uj E ) = —2— — - , (31) 

8 {uj e + Vpq^yi- 

where n a/ 3(q, oj e ) is the polarization function of the 
gauge field in Euclidean spacetime, and v E = xJa/-\/2ft 



is the Fermi velocity at the Dirac node. Hence, 



r(gaugc) 



oc — q Im 

A 2 g 1 



(v F q 2 — Auj 2 + irj) 1 / 2 



q 2 0(Au - v F q) 
(Aw 2 - V 2 F q 2 ) 1 / 2 



(32) 



For herbertsmithite, v F is estimated 1 ^ to be 5.0 x 
10 3 m/s. Hence, even at back scattering and with op- 
tical light at wavelength A « 500 nm, v F q corresponds to 
a frequency shift of approximately 4.0 cm -1 only, which 
is too small to be resolved by current instruments. There- 
fore, the collective excitation associated with the gauge 
flux will appear as a characteristic 1 /uj singularity in ex- 
periments. 

Physically, if the gauge boson is non-dissipative, the 
Green's function would have a simple pole, corresponding 
to a sharp delta- function- like signal. That we have al/w 
singularity in place of a delta function tells us that the 
gauge photon mode is strongly dissipative and is in fact 
overdamped. Following the analogy with the ordinary 
plasmon mode as stated above, this dissipative behavior 
of the emergent gauge boson can be thought of as the 
analog of the Landau damping. 
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VII. DISCUSSION 

In the previous sections we presented calculations of 
Raman intensity profile based on the Shastry-Shraiman 
formalism, assuming the validity of the U(l) Dirac spin- 
liquid state. We found a broad continuum in the Raman 
intensity profile in all symmetry channels, each displays 
a power-law behavior at low energy. Moreover, the pro- 
files are found to be invariant under arbitrary rotations 
in the kagome plane. For the E g and the Ai g channels, 
the continuum is accompanied by occasional sharp peaks 
that can be attributed to the various features of the DSL 
bands. In addition, the Raman intensity profile in the 
A-ig channel also contains a characteristic 1/cj singular- 
ity, which arose in our model from an excitation of the 
emergent U(l) gauge field. 

However, several caveats in our theoretical predictions 
should be noted. First, in order to compare with experi- 
mental data, \J must be converted to physical units. In 
Ref. [lH Hermele et al. estimated \ by fitting the spec- 
trum of projected one-particle excitations to the mean- 
field band structure, and found that x ~ 0.40. Together 
with J w 190A', this gives \J ~ 56 cm -1 . However, 
there is considerable uncertainty in this estimation, and 
so it may be a good idea to take x as a fitting parameter 
when comparisons with experiments are made. 

Second, when the contribution to the Raman inten- 
sity profile by spinon-antispinon pairs are calculated, the 
excited spinons and antispinons have essentially been 
treated as free fermions. However, they should really be 
regarded as complicated composite fermions, which inter- 
act with each other through an effective gauge field. Con- 
sequently, the actual excitation spectrum of the quasi- 
particles will almost certainly look quantitatively differ- 
ent from the ones presented here. Specifically, there may 
be finite lifetime effects, particularly prominent at high- 
energy, that causes the Raman intensity profile to be 
"washed out" compared with the ones presented here. 
Because of this, the sharp peaks that appear in Figs. [5] 
and [T^] may not be present in the actual data. Further- 
more, while the low-energy power-law structure of the 
Raman intensity profile is expected to survive, the de- 
tailed exponent is almost certainly modified from their 
mean-field values. Similarly, the contribution of Raman 
intensity in the A2 g channel by the emergent gauge bo- 
son may scale as J^ gauge ) cx (Aw) a , with a modified from 
-1. 

Third, in our derivation of the operators that corre- 
spond to the Raman transitions in the different channels, 
Eqs. [T0TU3| we stopped at the zeroth order in t/U and 
the leading order in t/(u>i — U). While terms higher or- 
der in t/U can be safely neglected, the same cannot be 
said for tj {uji — U), particularly near resonance (uj ~ U). 
Therefore, these higher-order contributions, which mod- 
ify the Raman intensity profile from those presented in 
the previous sections, may show up in actual data. In 
particular, the Raman intensity profile may not exhibit 
the abrupt drops as in Figs.lBland llQI Furthermore, since 



-pair c< Aw, a power-law with exponent much closer to 
1 may be found in the E g and Ai g channels at low energy. 

Fourth, while we argued that al/w singularity should 
be present in the A 2g channel, we have lost track of the 
ratio between its contribution and that by the spinon- 
antispinon pairs. Since the intensity of this singularity is 
proportional |q| 2 , it may be difficult to detect in optical 
Raman spectroscopy, in which the momentum transfer q 
is much smaller than the inverse of lattice constant. 

In the introduction, we mentioned the VBS state as a 
close competitor of the DSL state, as well as the proposal 
by Cepas et al. to distinguish between the two using the 
angular dependence of the Raman intensity profile. In 
light of the results presented in this paper, there is an- 
other feature in the Raman spectrum that can be used 
to distinguish between the two states. For in general, 
in a VBS state the spin excitations is gapped, while in a 
spin-liquid state it is gapless. Consequently, the Raman 
intensity profile should show an exponential dependence 
in the former case, and a power-law dependence in the 
latter case. 

Recently, Wulferding and Lemmcns^ have obtained 
Raman intensity data for herbertsmithitc. Their data, 
extending from 30 to 1500 cm -1 , shows a broad back- 
ground that persists beyond 500 cm -1 , in addition to 
a quasielastic line and several sharp peaks at finite fre- 
quency shifts. Furthermore, at low temperature (5 K) the 
quasielastic line is suppressed and the low-energy portion 
of their data shows a linear dependence with respect to 
the Raman shift. While their data are quantitatively 
different from the results of our theoretical calculations 
presented in Sees. IIVI IVI1 the existence of a broad con- 
tinuum can be seen as consistent with the U(l) Dirac 
spin-liquid model, even though the appearances of the 
other features would require the consideration of extra 
contributions (e.g., from the Zn impurities^^) that are 
not present in our model. 
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Appendix A: Derivation of the Raman transition in 
the Ai g channel 

In this section we shall consider the derivation of Ra- 
man transition rate in the e^ef — e^ef channel in more 
details. For completeness, we shall present derivations 
not only for the kagome lattice, but also for the square, 
the triangular, and the hexagonal ones. Although the 
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irreducible representation that corresponds to the polar- 
ization e^e^ — £/ e f may be named differently in these 
lattices, we shall abuse notation and continue to refer to 
them as the "A 2g " channel. 

We shall show that, contrary to claim by Shastry and 
Shraiman, the matrix elements in the A 2g channel van- 
ishes up to the t 4 1 (uji — U) 3 order in the square lattice 
according to their formalism. 

To extract the A 2g channel from the general polariza- 
tion matrix, note that given any particular hopping path- 
way, a "reversed pathway" can be constructed, in which 
all electron operators are conjugated and their order re- 
versed [for example, {c\c2){c\c^)(c\ci) is the reversed 
pathway of (c{c 3 )(cgC2)(c2Ci)]. Then, ejef i-> e^ef and 
the order the spin operators thus obtained are inverted. 
Hence, to the i 4 /(wj — U) 3 order, which corresponds to 
at most four spin operators, the only terms that sur- 
vive in the A 2g channel arc the spin-chirality operators 
Si ■ (Sj x Sfe). Thus, in our derivation it suffices to extract 
the spin-chirality contributions from pathways whose ini- 
tial and final currents are not co-linear. 

To depict the hopping pathways efficiently, the follow- 
ing abbreviations are introduced in the diagrams. A thick 
(blue) arrow is used to indicate the initial or the final hop 
in which a holon-doublon pair is created or destroyed. 
For the internal hops, the movement of a doublon is indi- 
cated by a thin (magenta) unbroken arrow and the move- 
ment of a holon is indicated by a thin (magenta) broken 
arrow. Lower case roman letters are used to indicate 
the ordering of hops. Note that in this scheme, a solid 
magenta arrow from i to j corresponds to the electron 
operators (ejej), while a broken magenta arrow from i to 

j corresponds to the electron operators (cjcj). 

The lowest order at which the spin-chirality term can 
show up is t 3 /(oji — U) 2 , which corresponds to pathways 
with one internal hop. Such pathway can be found in the 
triangular or the kagome lattice, or when next-nearest 
hopping is included. We shall show that the contributions 
to the A 2g channel by these pathways cancel in pairs at 
this order. 





(a) 



(b) 



FIG. 13: Two types of one-internal-hop pathways. Thick 
(blue) arrows denote initial or final hops in which a holon- 
doublon pair is created or destroyed, thin (magenta) unbroken 
arrows denote the movement of doublons, and thin (magenta) 
broken arrows denote the movement of holons. Lower case 
roman letters are used to indicate the order of hops. In (a) 
the internal hop is performed by the doublon while in (b) the 
internal hop is performed by the holon. 



It is easy to see that there are two types one-internal- 



hop pathways in general, both involving three lattice 
sites. In a pathway of the first type, a holon-doublon pair 
is created across a bond by an incident photon. Then, 
the doublon moves to a third site before rccombining with 
the holon to emit a Raman-shifted photon [Fig. 13(a)| . 
A pathway of the second type is similar, except that it is 
the holon that moves to a third site before recombining 
[Fig. [13(b)]. 

Applying the procedures as explained in Sec. ITT1 the 
operator that corresponds to the pathway in Fig 
is given by: 

(e/ • xi 3 ) (ej • x 2 i ) j 



13(a) 



x (c t ic 3 )(4c 2 )(c^ci) 

^13^32^21 



( e / ' x i3)(e 
(e/ -xi 3 )(e 



x 2 i) 
x 2 i) 



(w< - U f 

^13^32^21 



tr{X3X2Xi} 



2i S 3 • S 2 x Si 



where x 



(Al) 

x, — Xj is the vector from site j to site i, 
and "=" denotes equality upon neglecting terms that do 
not contribute to the A 2g channel. 

Similarly, the operator that corresponds to the path- 
way in Fig. |13(b)| is given by: 

T\,h = (e/ ■ x 3 i)(e. t ■ X12)- 



(oji - uy 



(4ci)(4c 3 )(clc 2 ) 



(e/ ■ x 3 i)(e 4 ■ X12) 
(e/ • x 3 i)(e 4 • x 12 ) 

—Ti e , 



^ 3 1^2 3 il2 

fa - uy 

£ 3 ii2 3 ii2 
(Wi - Uf 



(-l)tr{xiX2X 3 } 



2iSi • S 2 x S 3 



(A2) 

where Uj are assumed to be real in the last step. 

Since the two pathways depicted in Fig. [13] always 
come in pair, the contribution to the A 2g channel by 
one-internal-hop pathways vanishes upon summing as 
claimed. 

Now consider pathways that involve two internal hops, 
starting with the square lattice. Henceforth we shall as- 
sume that hopping is between nearest neighbors only, 
uniform, and real. The abbreviations C2 = t 4 /(uji — U) 3 
and Si-j-k = Si ■ (Sj x Sfc) will also be used. 

To count the two-intcrnal-hop pathways in the square 
lattice systemically, we fix the initial holon at site 1 and 
the initial doublon at site 2, and align the coordinates 
so that y — for site 1 and 2 and that X21 = x. All 
other pathways are clearly related to the ones satisfying 
the above conditions via symmetries. For the final hop 
and the initial hop to be non-collinear, a third site not 
collinear with site 1 and 2 must be involved, and we may 
further restrict our attention to pathways in which the 
third site lies in the y > half-plane, since the remaining 
pathways are related to these via the mirror reflection 

y -> -y- 
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There are four pathways that satisfy the above restric- 
tions, which are precisely the ones depicted in in Fig. [1] 
Applying the procedures as explained Sec. HH the contri- 
butions by these pathways are given by: 



T 2 , a = C 2 ef (-ep(c[c 3 )(4c4)(4c 2 )(4ci) 
-f 



= -C 2 efe u f tr{x3X4X2Xi} 
= -jC 2 e^e^(53 ; 4;i + 5 3;2 ;i + 5 4;2; i - 5 3; 4 ;2 ) , 



T 2 .b = C 2 e^^( C ic 2 )(4c4)(c t lC3 )(c t 2 c 1 ) 



C 2 e^e^tr{x 2 XiX3X4} 
zC 2 ef ey(«Si;3;4 — 5 2;3; 4 — 6>4 ;2; i — iS 2 ,i ; 3) , 



(A3) 
(A4) 



To 



C 2 e^(-e-p(4 C 4)( C t lC2 )(4c 2 )(4c 1 ) 
-C 2 efe«((c 4 4)(4c 2 )(4c 1 ) + (c 4 4)(4c 2 )( 

( tl '{(- 1 )X4X 2 Xi} + tr{x4}tr{x 2 Xi} 



eld 



= -C 2 efe^ 
= -2iC 2 ie*e v f S 4 ;2;i 



-i 



(A5) 



T 2 , d = C 2 e?e^(4 Cl )(4 C2 )(4c3)(4ci) 

= C 2 ere^((4c 2 )(4c 3 )(4ci) + (c\c,){c\c 2 ){c\ Cl ) 

= C 2 e*e) ( tr{ (-l)x3X2Xi} + tr{x 3 } trfexi}) 
" 2iC 2 ie* e^ 3;2; i . 

(A6) 

Summing all four terms, we found that T 2 a + T 2 ^ + 
T 2 . c + T 2 ,d = 0. Hence, for the square lattice with only 
nearest-neighbor hopping, the operator that corresponds 
to the Raman transition in the A 2g channel, Oa 2s > van- 
ishes to the t 4 /(uji — U) 3 order. 















y 






























1 



x — > X 

y — > —w 




FIG. 14: Mapping between the pathways in the square lattice 
and the pathways in the triangular lattice. 



Notice that the orthogonality between the x and y has 
not been invoked in the above derivation. Consequently, 
the above derivation carries to the triangular lattice upon 
mapping x and y in the square lattice to any two of bond 
directions in the triangular lattice. See Fig. [T3] for il- 
lustration. It can be checked that all two-internal-hop 
pathways with non-collinear initial and final hops in the 



3 -aa 



J — ai 




FIG. 15: The honeycomb lattice (thin gray lines), wherein the 
site label i, j, the unit vectors x, v, w (thin black arrows), and 
the primitive lattice vector (thick black arrows), are defined. 




J - 32 




(a) 



(b) 




% + a 2 




3 ~ ai 



(c) 



(d) 



FIG. 16: Pathways with two internal hops that contribute 
to the A-zg channel in the honeycomb lattice, with the initial 
holon fixed at site i and initial doublon fixed at site j. 



triangular lattice can be obtained from such mappings 
and that there is no issue of double-counting. Hence, we 
conclude that Oa 2s vanishes up to the t 4 /(a;, — U) 3 order 
in the triangular lattice also. 

Evidently, the criterion that any three non-collinear 
nearest-neighbor sites belong to a four-site loop is a cru- 
cial ingredient for the cancellations of the two-internal- 
hop contributions as seen above. This criterion is not met 
in the honeycomb lattice or in the kagome lattice. Hence, 
Oa 2s ma y not vanish in these lattices at the t 4 /( w i — U) 3 
order. 

First consider the honeycomb lattice, which is shown in 
Fig- El wherein the site labels i, j, the unit vectors x, v, 
w along bond directions, and the primitive lattice vectors 
ai, a 2 , are defined. For an initial holon at i and an initial 
doublon at j, there are four two- internal-hop pathways, 
listed in Fig. \Wi Summing up their contributions, we 
get: 



T 2 ,i.j = 2iC 2 (^ — ef ej<S,+ ai ;j;i 



(A7) 
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If the initial doublon is fixed at j — a 2 or j — ai instead, 
the contributions arc, respectively, 



^2,i,j— a2 — 2iC*2 ^ 6i&fSi—a2;j~EL2Vi ~b ^ f ^j— a i;i _ a 2 ;i 

(A8) 

2l2,ij'-ai = 2iC 2 ^ — e^eJiSj+ a2 _ ai; j_ ai; j + e"ej;5j ;:) _ ai; j 
— e^Oi_ ai;:) _ ai; i e i efOj- a2 -j- ai ;i J ■ 

(A9) 

And the analog of Eqs. lA7l IA91 when the holon is fixed 
at j, are given by: 



1 2 ,., = 2iC 2 ( - efe'?S j - ai . j . i + e?e}S, 



«+a 2 ;j';i 



' *,'< • a y:- ) , 



(A10) 



2 ,j,i+a 2 



2iC 2 \ e i ey<>j-|- a2; .;+ a2; j + e,j i5j-|_ ai; j_|_ a2; j 

e i t 5j+a2-ai;i+a 2 ;j T" c f<->i:i+a. 2 ;j 



(All) 

n2,j,i+ai = 2iC*2^ — e ™ eJ5j+ ai — a2 ;j-|- ai; j + ej"e^5i ; ;-|- ai ; j 



°i °/'-'j+ai ;i+ai ; 



-i e / i+a 2 ;i+ai;j 



(A12) 

Summing Eqs. IA7IIA121 over all lattice vectors {R}, 
and reorganize slightly, we finally obtain: 

/2 ,„x = UC 2 ((£,■*< a + S i;i;i+ai )(e?e? - ej"e?) 

R 

+ (<5i-a 2 ;i;j + <S i+ a 2 « - e?g}) 

+ (^-a lii; i-a 2 +5 i+ai;i;i+ a 2 )(efe} - e?e?)) 

= ^ e (( + (^7 - *m 



(A13) 

where graphical symbols are introduced on the second 
equality to denote the spin-chirality operators. Note that 
even though the site labels are omitted in the symbols, 
upon the summation over lattice vectors {R} there is no 
ambiguity as to which spin-chirality operator a particular 
symbol is referring to. 

1 x _|_ VjL 
2 C T 2 



Using e v 



Eq. IA13I can be converted back to the Cartesian coordi- 




I + aj. — a 2 



FIG. 17: The kagome lattice (thin gray lines), wherein the site 
label i, j, the unit vectors x, v, w (thin black arrows), and 
the primitive lattice vector (thick black arrows), are defined. 





(a) 



(b) 



- ai. - a 2 




(c) 



(d) 



FIG. 18: Pathways with two internal hops that contribute to 
the A 2g channel in the kagome lattice, with the initial holon 
fixed at site i and initial doublon fixed at site j. 



nates, which yields: 



+ 



(A14) 



Finally consider the kagome lattice, which is shown in 
Fig.[17l wherein the site labels i, j, £, unit vectors x, v, w, 
and primitive lattice vectors ai, a 2 , are defined. Fixing 
the initial holon at site i and the initial doublon at site 
j, four two-internal- hop pathways contribute to the A 2g 
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channel (Fig. [T8|) . The sum of their contributions is: 



T, 



2,i,j 



„x-v c i a x z w C 

~ i /°^+ai — a 2 ;j;i e i e / a 2 ;j; 



(A15) 




where graphical symbols are again introduced on the sec- 



J 



ond equality. Note again that upon summing over all 
lattice vectors {R} there is no ambiguity as to which 
spin-chirality operator a particular symbol refers to. 



By changing the site where the initial doublon is lo- 
cated, we get, upon summation, the following contribu- 
tion to Oa 2 by two-internal-hop pathways whose initial 
holon is located at site i: 



T 2 j = 2iC 2 [e?e) 



^(A+V 



■;(A+V 



+ 




(A16) 



Obtaining the contributions to Oa 2s by two-internal-hop pathways whose initial holon is located at site j or £ in 
an analogous manner, we finally get, upon summing over all lattice vectors and basis sites, 



T 2Mg = 4*C 2 ]T ( (ef ej - efe* f ) 
R V 



Or, converting back to the Cartesian coordinates: 
T 2Mg = 2V3iC 2 ^>?e? ~ e?ty ( 3 A +3 V 

R 

which is what we quoted in Eq. 1131 




(A17) 




(A18) 



In summary, we found that, with only nearest- neighbor 
hopping, the A 2g channel Raman T-matrix does not van- 
ish in the honeycomb lattice or the kagome lattice at the 
t 4 /(U — uji) 3 order, but does so to this order in the square 
lattice and the triangular lattice. 



Appendix B: Derivation of Oa 13 and Oe s to the 
t 4 /(LUi - Uf order. 

In this section we shall derive Oa 19 and Oe 9 to 
the t 4 / (oji — U) 3 order for the kagome lattice. As al- 
ready noted in Sec. HH at the t 2 /(uji - U) order the 
Shastry-Shraiman formulation reproduces the Fleury- 
London Hamiltonian. In the A\ g channel this gives rise 
to an operator proportional to the Hciscnberg Hamilto- 
nian, and in the E g channel it gives rise to the operators 



O (i) and O (2) as shown in Eqs. [TO1[TT1 

At the t 3 jiiOi — U) 2 order, it can be checked that the 
two types of pathways depicted in Fig. [T3] cancel each 
other not only in the A 2g channel but also in the E g 
and A\ g channels. Thus, it remains to consider path- 
ways having two internal hops, which contributes at the 
t 4 /(w ?; - U) 3 order. 

By considering pathways and their "reserved" coun- 
terparts, in which all electron operators are conjugated 
and their order inverted, it is shown in Appendix [X] that 
only spin-chirality operators contribute to the A 2g chan- 
nel. From the same construction, it can be seen that the 
spin-chirality operators do not contribute to the E g and 
A\ g channels. Moreover, we are interested in the inelas- 
tic, and hence non-constant, part of the Raman transition 
operators. Therefore, to determine the t A / (oji — U) 3 or- 
der terms in Oe b and O A 2g > it suffices to extract the spin 
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(e) (f) (g) (h) 



FIG. 19: Pathways with two internal hops in a kagome lattice, with the initial holon fixed at site i and initial doublon fixed at 
site j. 



dot product terms for each process. 

Fixing the initial holon at site i and the initial doublon 
at site j, there are eight more two- internal- hop pathways 
that contribute to the Ai g and E g channels in addition to 
the four depicted in Fig. [THJ These are listed in Fig. [121 

Applying the procedures as explained in Sec. [TH the 
spin dot products resulting from Figs. 18(a)] [l8(b)[|19(c~ 
and |19(d"J| are, respectively: 
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kg.a = C 2 ef(-e v f )(^tr{(-l)xeXjXt} + tr{xe}tT{x 3 Xi} 
= C 2 eU~e})(S e -Si -Si- Sj - Sj • S/) , 

(Bl) 



2k, 



kag,6 = C 2 eU-ey)[tr{(-l)xiXjXi}+Hxe}tr{ Xj Xi} 



C 2 eU-ef){S 3 -St -Si- Sj - S t • S<) , 



(B2) 



T kag ,c = C 2 eU-e x f )(4cj)(clc e )( C l Cj )(c]c t 
= C 2 ef(-e})tT{ Xe }tr{x jXi } 



2C 2 er(- 



zJ)Si ■ Sj , 



= C 2 ef(-e^)tr{x^}tr{x J -Xi} 
= - 2C 2 e?(-epS i • S, , 



(B3) 



(B4) 



where = denotes equality upon neglecting additive con- 
stants and spin-chirality terms. The contributions of 
other pathways in Figs. [18] and [19] can be obtained by 
relabeling sites and vectors that appear in Eqs. IB1IIB4I 

The sum over all pathways in Figs. [TBI and [191 gives: 



Tkag,i,j — C 2 ^ e f e /(Sj: • Sj_ ai + Sj ■ Si+ ai + 14S; • Sj Si ■ Si+ai Sj ■ Sj_ ai ) 

+ e*e V f{Si ■ S£+ ai -a 2 + Sj • Si — Si ■ St — Sj ■ Sf +ai _ a2 ) + efej(Sj ■ S^_ a2 + Si ■ Si — Sj ■ Sg — Si ■ S<>_ a2 )^ 

(B5) 

The contributions when the initial holon and/or the initial doublon are located at other sites can be obtained by 
relabeling. The sum over the locations of the initial holon and the initial doublon yields: 



Hag = 2C 2 J2 {' >' /( l,iS < ■ Sj + 16Sj ■ S i+ai - 2Si ■ S i+ai - 2Sj ■ Sj_ ai ) 
{R} V 

+ e^ej(l6Sj • Si + 16Sg ■ Sj +a2 _ ai — 2Sj • Sj +a2 _ ai — 2S^ • S^ +a2 _ ai ) 

+ e™e™(l6Si • St + 16Si ■ S i+SL2 - 2S t ■ S i+SL2 - 2S e ■ S e+a2 ) 

+ ( e i&f + e r^/)(Sj ■ S£_ a2 + S( ■ Sj_ ai — Sj • S( — S^_ a2 • Sj_ ai ) 

+ ( e f^/ + e "^)(Si ■ S(> +ai _ a2 + S( ■ Sj +ai — Si ■ S( — S£ +ai _ a2 • Sj +ai ) 

+ ( e i^J + e r6/)(Si • Sj +a2 _ ai + Sj • Si +a2 — S; • Sj — Sj +a2 „ ai • S 1+a2 )^ 



(B6) 
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Projecting onto the A\ s channel and neglecting a piece proportional to the Heisenberg Hamiltonian yields Eq. [12l 
while projecting into the E„ channels yields the following contributions to 0„d) and 0„m at the t A /(uii — U) 3 order, 

&g i^g 

respectively: 



E {i) = C2^(2S; ■ Sj+a 2 -ai + 2Sj ■ S; +a2 — Sj ■ Sf_ a2 — • Sj_ ai — Sj • S£ +ai _ a2 — Sg ■ Sj_|_ ai ) 

(2Sj • Si-)_ ai + 2Sj • Sj_ ai Si ■ Sj-|_ a2 • S^-|- a2 Sj ■ Sj-|_ a2 _ ai S^ • S£+ a2 _ ai ) 
+ 7(2Si • Sj + 2Sj ■ S i+ai — Si -S? — Si ■ S i+EL2 — Sj ■ Si — S? ■ Sj +a2 _ ai )^ , (B7) 

80 E (2) = V3C2^{Sj ■ S£_ a2 + S( ■ S 3 _ ai — Si ■ S(> +ai _ a2 — S( ■ Sj +ai ) 

— (Si • Si+ a2 + Sg ■ Sf+ a2 — Sj ■ Sj+ a2 _ ai — Si ■ Sf+ a2 _ ai ) 

+ 7(Si ■ Se + Si ■ Si +a2 — Sj ■ Se — Se ■ Sj +a2 _ ai ) ) . (B8) 



I 

As explained in Sec. IIV1 it can be checked that SO„(i> transition in the Dirac Hamiltonian Eq. [20l 
and 80„(2) in the above equations do not produce any 
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